The fluctuation-driven transition between metastable states is particularly rel- 
evant to many important events in physics, chemistry, biology, etc. Examples 
include chemical reactions, biological switches, nucleation processes, to list just 
a few. In many systems, the dynamics may involve delayed interactions due 
to limit transmission rate of matter, energy or information transport, or some 
kinds of feedback. Thus, it is of great interest to investigate the impact of de- 
layed interactions on the transition behaviors of such systems. In this paper, 
we address this issues by investigation the effect of delayed interaction on the 
dominant transition pathways which is informative to explain the mechanism of 

C/3 . fluctuation-driven transitions. 
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O ! We investigate delay effects on dominant transition pathways (DTP) between 



metastable states of stochastic systems. A modified version of the Maier-Stein model 
with linear delayed feedback is considered as an example. By a stability analysis of 
the 'on-axis' DTP in trajectory space, we find that a bifurcation of DTPs will be in- 
duced when time delay r is large enough. This finding is soon verified by numerically 
derived DTPs which are calculated by employing a recently developed minimum ac- 
tion method extended to delayed stochastic systems. Further simulation shows that, 
the delay-induced bifurcation of DTPs also results in a nontrivial dependence of the 
transition rate constant on the delay time. Finally, the bifurcation diagram is given 
on the t — /3 plane, where /3 measures the non-conservation of the original Maier-Stein 
model. 
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I. INTRODUCTION 



Real dynamical systems are often subject to weak random perturbations, such as thermal 
noise at a nonzero temperature. It has been a common sense that small fluctuations can 
produce a profound effect on the long time dynamics by inducing rare but important events. 
For instance, fluctuations may result in transitions between metastable sets of determinis- 
tic dynamical system, which can be related to a large number of interesting phenomena in 
physics, chemistry and biology such as nucleation processes, chemical reactions, and biolog- 
ical switches. 

In recent years, fluctuation-driven transitions (FDT) have gained great research attentions^-. 
One of the fundamental purposes of studying FDT is to explain how the transition occurs. 
Freidlin-Wentzell theory of large deviations provides one of the right frameworks to under- 
stand FDT 10 — 2 . When the amplitude of fluctuation is small, the distribution of trajectories 
which make transitions between metastable sets is often sharply peaked around a certain 
deterministic path or a set of paths. It then becomes very important to identify such 
dominant transition pathways (DTP) which can be highly informative to help elucidate the 
underlying mechanism of the FDT. Usually, the DTP tells how the transition happens step 
by step, identifies the transition state(s), and can also be used to derive other important 
quantities such as the transition rate of the FDT. In conservative systems, where there 
exists an underlying energy landscapes^!, the DTP is actually the minimum energy path 
and is everywhere tangent to the potential force. In this case, the DTP must first approach 
a transition state which is usually a saddle point on the basin boundary of one attractor, 
and then runs along the unstable manifold of that point and enters the basin of attraction 
of the second attractor. The DTP before reaching the transition state is actually the time 
reversed heteroclinic orbit of the unperturbed system joining the attractor and the saddle 
point. In non-conservative systems in which detailed balance is absent, some new interest- 
ing phenomena have been found, e.g., symmetry breaking bifurcation of the optimal escape 
path can be observed", an unstable fixed point 2 or an unstable limit cycled can act as the 
transition state. Moreover, a complex transition paradigm containing a saddle point and 
two limit cycles as transition states was reported in the Lorenz system^. The DTP has also 
been used to explore the configuration space of systems with complicated structure^, and 
study the nucleation process in the presence of shear in a two-dimensional Ginzburg-Landau 
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equation 4 . 

Nevertheless, most studies of DTP so far are limited at least in one sense, i.e., the state 
of system at one time can only influence and be influenced by its state at that same time. In 
fact, a variety of sources, such as the limit transmission rate of matter, energy or information 
transport, or some kinds of feedback, might allow events at one time to affect the state of 
the system at some later time. In these cases, time delayed variables and equations should 
be used to describe the dynamics. As we already know, delay models have been widely used 
to describe chemical kinetics 1 ^, neuronal networks 1 ^, circadian oscillators^"—, physiological 
systems^, optical devices^, and so on. In addition, time delay can lead to a variety of 
interesting and important phenomena, such as delay-induced oscillation 2 ^, delay-induced 
excitability 21 , delay- induced oscillator death 22 . However, to the best of our knowledge, 
how delay would influence the FDT dynamics of a stochastic system, albeit its apparent 
importance, has not been studied yet. 

In present paper, we have addressed such an issue by investigating the effect of time delay 
on the DTP in a modified version of the Maier-Stein model with linear delayed feedback. By 
an analysis using small delay approximation, we find that the DTP undergoes a bifurcation 
via transverse instability in the trajectory space when the delay time r bypass a certain 
threshold value r c . By extending a recently developed minimum action method 23 to this 
delayed stochastic system, we have also obtained the DTPs by numerical simulations, which 
further confirms the analytical results. In addition, this bifurcation of DTP results in a 
nontrivial phenomenon of the FDT: The transition rate constant between two metastable 
states shows distinct dependence on the delay time below and above the bifurcation point. 
Finally, the bifurcation diagram is given on the r — j3 plane, where (3 stands for the non- 
conservation effect of the original Maier-Stein model. 

II. ANALYSIS 

In general, a delayed stochastic system whose dynamics is determined by both the present 
state x(t) and the state x(t — r) with the delay time r > can be described as 

x(t) = F(x(t), x(t - r)) + y/£a(x)ri(t), (1) 

where F(x(t), x(£ — r)) is a known drift vector field and r](t) is a set of independent Guassian 
white noises with zero mean and unit variance, e is a small positive number, and cr(x) is 
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FIG. 1. The potential field and the minimum energy path (the bold solid line) of the original 
Maier-Stein model. 

related to the diffusion tensor by a = oo T . 

To show the effect of delay, here we consider a modified version of the Maier-Stein model 
with state vector x = (u, v) as an example, whose linear term is modified to be a delayed 
feedback. 

ii(t) = ut-r -u 3 - (3uv 2 + y/er] u (t) 

(2) 

v(t) = -v t -T - u 2 v + s/er] v (t) 

When r = 0, Eq.(|2]) recovers to the original Maier-Stein model. In the absence of noise 
terms, it has two stable steady states at A = (—1,0), B = (1,0) and a saddle point at the 
origin (0, 0) for all values of /3 > 0. In the presence of weak noise, however, both A and B 
become metastable states (MSS), and the FDT from one MSS to another is a rare event. 
reflects the non- conservation of the original model. For /3 = 1, the drift field of the original 
system can be viewed as a gradient of a potential field, and the DTP from A to B (or vice 
versa) is actually the minimum energy path connecting A and B along the u axis, which is 
shown in FigHJ 

In the delayed Maier-Stein model whose dynamical equation is governed by Eq.([2J), A and 
B are still the asymptotic fixed points of the system. By a simple linear stability analysis, 
we can find that A and B are stable for r > 0. Similar to the case without delay, when small 
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perturbation is present, FDT from one MSS to another is allowed. As this study focuses 
on the effect of delay, we will fix = 1 if not otherwise stated to avoid the influence of 
instantaneous non-conservative effects. 

We now try to figure out the critical value of r at which the on-axis DTP will be unstable. 
To this end, we expand Eq.flSJ) in powers of r using a small delay expansion around x(i)2£ 
as follows, 

u(t) = (1 - r)(w - w 3 - uv 2 ) + y/e(l - r)r] u (t) 

(3) 

v(t) = (1 + t)(-v - u 2 v) + v^(l + r)r] v (t) 

When weak noise presents, for a given transition path $ = {$ = A, ...,$t = B} from 
MSS A to B over a finite time interval T, a Freidlin-Wentzell action functional Sy[$] can 
be calculated by path intergral along $. The main result of the Freidlin-Wentzell theory^ is 
that for sufficiently small e, a probability can be assigned for path <f> as P(|{x} — $| < e) w 
exp{—S T [^]/e}, where e is a sufficiently small positive number. The probability suggests 
that DTP(s) is(are) the transition path(s) minimizing the action functional. For the on-axis 
DTP ip = {tp = A, ...,ipT — B,ip v = 0}, the Freidlin-Wentzell action functional St[iP] can 
be written as 

\q-G\ 2 dt , (4) 
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where F a = ((1 -t)(u-u 3 -uv 2 ), (l + r)(-v -u 2 v)) = (F a u , F%), q = (u/(1-t),v/(1 + t)), 
G(u,v) = (F a u /(l-r),F a 7(l + r)) and r ^ 1. Here, we use |q| 2 + |G| 2 > 2|q||G| to get the 
inequality. Note that, this inequality is actually an equality, since the DTP ip minimizes the 
action functional. The W(u,v) is given by 



-(u,v) 

W{u,v) = I -G-qdt 



' A i ' „,2 , 1 22 



(5) 



2(l-r) v 2 ' ' 2(1 + r)" ' 1 - r 2 

The right-hand side of the inequality in Eq.(jl} is a line integral along the directed curve 
tp, which can be considered as a geometric action functional S similar to the one in Ref.-. 



S = Jq |q||G|d£ = • dq, where 9 is the angle between G and q. As both G and q 

are always along the u axis, cosO = 1 or — 1 for any point (u,v) G ip. Then, we can expand 
the geometric action functional of the segment from A to (u, v) near the u axis in powers of 
v as 



S[{^o,...,ip t = (u,v)}} 



\W(u,v) - W{A)\,u< 
\W(u,v)\ + \W(A)\,u> , (6) 
= m (u) + m 2 {u)v 2 + o(v 2 ) 

where rri2{u) is a measure of the transverse stability of the transition path at the point 
{u,v). If rri2{u) < 0, any small perturbation in the v or — v direction will lead to a smaller 
value of the action functional, and the DTP along the u axis will be unstable. Notice that 
the perturbations in both directions have no difference in decreasing the value of the action 
functional, which indicates that, there will be two equivalent DTPs, related by v — > —v, if 
m 2 {u) < 0. 

For the left segment of ip where u < 0, we have 

, . i 2(1+7) + (T^)^ 2 ' r < 1 /-x 
m 2 {u) = { y ' , (7) 

2(l+r) (1-r 2 ) " ' ' ^ 1 

Then, a straightforward calculation shows that the instability region of the on-axis DTP is 

r > 1. (8) 

Similarly, we find that, for any point with u > 0, m 2 (u) > 0, which means that small 
perturbations at the right segment do not affect the stability of the on-axis DTP no matter 
what value r takes. 

It is noted that, Eq.Q is derived under small delay assumption, and such an expansion 
has been shown to be valid to quadratic order in . As Eq. ([8]) is not small enough, we can 
not expect that the derived instability boundary is the exact one. Even so, Eq.([H|) suggests 
a guiding picture that, there will be a threshold above which the on-axis DTP will undergo 
a bifurcation via transverse instability on the left segment. What's more, as the geometric 
action functional 5* is independent on T, the bifurcation of DTP and instability condition 
Eq.()8]) will also be independent on T. 



III. NUMERICAL SIMULATION 



To verify the analytical result, we now derive the DTP by simulation using a recently 
developed minimum action method 23 . The extension of this method to a delayed system 
is straightforward except for some details. One should note that the delay time defines a 
upper limit of the step size At = T/N when we discretize the time domain [0, T] to a iV-size 
mesh. To consider the delay effect, At should be properly chosen so that At = r/m, where 
m > is a positive integer. We start from a test path $ with $o — A and $t = B, and 
update the path till convergence by iterating on solving the gradient dynamic as 

d<$> t 5S T [<$>] 



dk 8<$> t 




< t < T 

(9) 



dk dk 

Where, k > plays the role of pseudo time for the updating, and the action functional of 
the path St[§] can be calculated by the first equality of Eq.(@| by using $ and F instead 
of ip and F a respectively. The resulting path(s) with minimum action is(are) our DTP(s) 
ip, and the rate constant P that the transition from A to B occurs can be approximately 
estimated by 

P x lim exp{--SrP0]}, (10) 

where f(e) x g{e) if log f(e)/ log g(s) — > 1 as e — > 0. When the DTP bifurcates, P is 
calculated via summing ones of each path. 

The numerical DTPs for r = 0.0, 0.3 and 1.2, and T = 100 are plotted in Figj2j It can be 
observed that two equivalent DTPs, symmetric under v — > —v, appear, which means that, 
the DTP does undergo a bifurcation in the trajectory space when r is large enough. When 
t = 0.3, only one on-axis DTP exists. While r = 1.2, two off-axis DTPs are observed. For 
further quantitative analysis of the DTPs, we define the life time, as defined by X. Zhou, 
etc. 8 , to be the (signed) time in the deterministic system starting from a given position to 
reach the K-neighborhood of A (negative life time) or B (positive life time), where k = 10 -5 
in our simulation. A negative life time indicates that the point is in the attraction basin 
of A, and a positive one implies that the point will be attracted to B. Thus, the transition 
state iptran can be determined by the middle of the last point with negative life time and the 
first point with positive life time. The life times tuf e for each point on the numerical DTPs 




FIG. 2. The dominant transition pathways of delayed Maier-Stein model when r = 0.0, 0.3 and 
1.2. 



are presented in Figj3]Ja). We find that, all the DTPs are separated by the origin (0,0) into 
two segments which lie in basins of attraction of A and B, respectively. The fluctuation of 
the negative life time at r = 1.2 may due to the fact that, the trajectory back to A doesn't 
overlap with the DTP in off-axis region, and it is easier for some point to return to A in 
the deterministic delayed system. In addition, we have also calculated the amplitude of 



optimal fluctuation forces b, 



for each point on the numerical DTPs. 



optm = 1pt~ F(lp t ,lp t -T) 

The forces are considered to be optimal is due to the fact that, as they are calculated along 
the dominant transition pathway ip, we will get this dominant transition pathway back by 
applying these forces to Eq.(l). Since b optm is proportional to the deterministic force F 
which is at A and the transition state, its value will be near for the points close to A, 
then increases and reaches its maximum at some middle point, and decreases to nearly for 
the points close to the origin, as presented in FigB^b). Besides of this, several points can 
be addressed. Firstly, for all r, strong fluctuation force is needed for the system to escape 
from the attraction of A. Just after the system passes through the origin, b optm decreases to 
zero, immediately. Secondly, the curve at r = 0.3 is overlapped with its analog at r = 0.0, 
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FIG. 3. (a) The life time tuf e and (b) the optimal fluctuational force b op t m for each point of the 
dominant transition pathways shown in Figj2j 

which confirms that the on-axis DTP is still stable when the delay time is small. Finally, for 
r = 1.2, b op tm of the left segments show large discrepancies from the one at r = 0.0 obviously, 
but the right segments do not. The results given by Figf5]and Figj3]are consistent with our 
analytical results. 

In order to quantitatively describe the bifurcation of DTP, we introduce the maximal 
distance from DTP to u axis, L, as follows 

L = max(/ i ),0 < t < T (11) 

where It is the distance from ip t to the u axis. In FigJU^a), L as a function of delay time r is 
shown. It can be seen that L stays nearly zero for small delay time until r > r c w 1.1. We 
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FIG. 4. The bifurcation of the dominant transition pathway, (a) The maximal distance L from 
dominant transition pathway to the u axis as a function of delay time r. L arises from nearly zero 
to a remarkable value while delay time pass through the threshold r c ~ 1.1. (b) The scaling of 
threshold r c with observation time T. 



note that the threshold r c is close to the result given by Eq.([8]). To make sure the bifurcation 
is not a result of small T in our simulation, we test the scaling behavior of r c with T, which 
is shown in FigJU^b). The independence of r c on T implies that the bifurcation of DTP 
in trajectory space does occur for large Ts, which also confirms our analysis. For infinite 
T, it is not available to calculate DTP directly (A geometric minimum action method has 
been developed by Heymann et. al.— to deal with the infinite T problem, however, it is not 
suitable for delayed systems). 



IV. DISCUSSION 

An important quantity derived from the DTP is the action functional St[i/)] which can 
be calculated by the first equality of Eq.(j4j) by apply F instead of F a and is related to 
the transition rate constant by Eq. (fT0|) . St[iJ)] as a function of r is shown in Figj5](a). 
When t is small, St[iP] increases almost linearly as r increases, however, when r is large, 
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St[i/j] departs from the linear relationship obviously. To show this clearly, we plot the slope 
SSt[iP)/St in the insert of Figj5]^a). For r < 1.1, (5St[4>]/5t)/t w 0.0 which means the 
slope is nearly unchanged. When r > 1.1, the slope decreases as r increases at a rate of 
(SSt[iP]/St)/t « —4.5. For comparison, we also run the dynamic equation Eq.(J2J) directly 
involved with forward flux sampling approach (FFS)— to get the transition rate constant P 
at different r. P at e = 0.02 are plotted in Figfj^b). The nontrivial dependence of P on 
t is also observed, which is consistent with the result derived by DPT qualitatively. The 
quantitative difference is due to the fact that, while we estimate P by Eq. (|T0|) . the crossing 
through the transition state is considered as ballistic, i.e., it assumes that every crossing 
gives rise to a successful transition. For the diffusive crossing, some crossings may turn back 
to A and P is overestimated. The true P should include a prefactor Co, which evaluate 
the ratio of successful transition to total crossing, as well as the exponential factor Eq.f fTOj) . 
Since FFS simulation samples the transition rate constant directly by system's dynamics, 
we can calculate Co by the ratio between the two curves in FigJU^b). Fig^c) presents Co 
as a function of r. The abrupt increasing of Co near the bifurcation point r c implies that 
the system may undergo some sort of critical behavior as a result of delay induce DTP 
bifurcation. 

In the above simulation, we have fixed — 1 to avoid the influence of non- conservative 
effects other than time delay. It has been reported that, bifurcation of the DTP also occurs 
when /3 is varied^. To understand the dependence of the DTP on the both parameters r and 
/3, a two-dimensional bifurcation diagram in the r — /3 plane can be plotted. By extensive 
MAM simulations at different (3 and r, we plot this diagram in Fig. ((6]). It seems that, the 
bifurcation curve asymptotically approaches the line (3 = 0.5, i.e., the symmetry breaking 
bifurcation of DTP doesn't occur for any delay time when (3 < 0.5. We have tried to 
figure out the underlying mechanism theoretically, however, a potential-like quantity similar 
to Eq.([5]) is not available for the non-conservation case /3 ^ 1, which may need further 
research. 

V. CONCLUSION 

In summary, the dominant transition pathway between two metastable states of a delayed 
stochastic system is studied. To show the effect of delay, we consider a modified version of 
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FIG. 5. (a) The minimal action functional Sr[Y>], (b) transition rate constant P and (c) its prefactor 
Co as a function of delay time r. The insert panel of (a) is the slope 5St\^\/5t as a function of 
delay time r. Here, e = 0.02 for the forward flux sampling (FFS) simulation. 
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FIG. 6. The bifurcation diagram of dominant transition pathways in r — /3 plane. 
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the Maier-Stein model with linear delayed feedback as an example, whose original model is a 
typical symmetric 'double well' system, to study the noise induced transition. Our analysis 
by small delay approximation shows that, time delay will induce a new DTP via a bifurca- 
tion in trajectory space when the delay time passes through a threshold r c . By employing 
a recently developed minimum action method, we can calculate the DTP by minimizing 
the Freidlin-Wentzell action functional for transition path, numerically. The numerically 
derived DTP confirms that the bifurcation does occur for r > r c = 1.1. Other details of 
DTP bifurcation are also verified by numerical results, including bifurcation via transverse 
instability, bifurcation at left segment, arising of two equivalent DTPs after bifurcation, etc. 
From the DTP, the transition rate constant can be derived, which shows distinct depen- 
dence on the delay time below and above the threshold. This dependence is also observed 
by directly running the dynamic. The bifurcation diagram is also investigated. Since time 
delay is an important factor in many real systems, we believe that the present study can 
shed new light on understanding the mechanism of fluctuation-driven transitions in exper- 
imental studies and open more perspectives on the study of fluctuation-driven phenomena 
in non-conservation systems. 
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